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A  dominant  height  growth  model  and  a  site  index  model  were  developed  for  Calabrian  pine  ( Pinus  brutia 
Ten.)  in  central  Cyprus.  Data  from  64  stem  analysis  in  32  temporary  plots,  where  Calabrian  pine  was  the 
only  tree  species,  were  used  for  modeling.  The  plots  were  selected  randomly  in  proportion  to  two  site 
types.  Four  difference  equations  were  tested.  The  evaluation  criteria  included  qualitative  and  quantitative 
examinations  and  a  testing  with  split  data.  The  difference  equation  of  Korf  showed  the  best  results  for 
all  data.  An  analysis  of  the  height  growth  patterns  among  sites  -  as  these  were  defined  from  the  selected 
equation  -  was  made  in  order  to  study  the  behavior  of  different  site  index  curves.  Results  indicated  the 
validity  of  a  common  height  growth  model  for  the  two  sites.  In  spite  of  the  irregular  height  growth  pattern 
observed  in  Calabrian  pine,  the  model  obtained  allows  us  to  classify  and  compare  correctly  Calabrian  pine 
stands  growing  at  different  sites. 

©  2011  Elsevier  Ltd.  All  rights  reserved. 
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1.  Introduction 

The  species  Pinus  brutia  Ten.  (Calabrian  pine)  is  a  characteristic 
species  of  the  eastern  Mediterranean  [12,26,27,34,64,83].  P.  brutia 
is  widely  extended  in  Turkey  and  far  Eastern  Greece,  secondarily  in 
the  Crimea,  Caucasus  coast,  Azerbaijan,  Iraq,  Syria,  Lebanon,  Crete 
and  Cyprus  [54]. 

It  is  an  important  species  for  multi-purpose  forestry;  construc¬ 
tion,  industry,  carpentry,  firewood  and  pulp  [7,1 0,61  ].  Additionally, 
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its  post-fire  regeneration  ability  makes  it  a  unique  forest  species 
within  the  fragile  Mediterranean  ecosystems  [68]. 

P.  brutia  generally  occurs  in  the  Mediterranean-type  climate  of 
hot  and  dry  summers  and  mild  and  rainy  winters.  This  species 
occurs  most  abundantly  in  the  semi-arid  and  sub-humid  zones 
[64].  The  species  has  strict  rainfall  requirements,  being  absent  from 
the  arid  bioclimatic  zone  and  rare  in  the  lower  semi-arid  zone, 
but  widespread  in  the  humid  zone  [63].  We  find  Calabrian  pine  in 
zones  with  a  mean  annual  precipitation  between  400  and  2000  mm 
[36,45].  P.  brutia  is  well  adapted  to  the  Mediterranean-type  cli¬ 
mate  in  several  physiological  and  morphological  characteristics 
and  is  a  drought  resistant  species  [17,19,33,44,45,53].  The  adapt¬ 
abilities  of  its  provenances  to  drought  vary  [32],  but  P.  brutia 
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generally  achieves  its  optimum  growth  in  rainy  regions  of  more 
than  900-1000  mm  mean  annual  precipitation.  Its  typical  elevation 
range  is  between  0  and  1 500  m  above  sea  level  [45].  The  geograph¬ 
ical  isolation  among  populations  implies  the  existence  of  several 
races,  with  different  ecological  characteristics  and  behavior  [54]. 

Calabrian  pine  is  a  fast  growing  tree  species;  its  mean  annual 
increment  in  plantations  can  be  over  10m3/ha  (at  a  spacing  of 
4.5  m2 ;  2222  trees  per  hectare,  site  class  I,  site  index  20  m)  between 
the  ages  of  25  and  35  years  [35].  New  stands  of  P.  brutia  generally 
begin  after  fires,  although  it  can  also  naturally  regenerate  without 
fire  [26,27,46,47,48,76,77].  Several  adaptations  greatly  contribute 
to  its  post-fire  regeneration  [20,46].  Seedlings  develop  rapidly 
growing  tap  roots  while  the  stems  grow  comparatively  slowly  in 
height,  so  that  seedlings  can  have  65  cm  tap  roots  within  five  to  six 
months  after  germination  [26,27,46].  Early  stem  growth  is  greater 
at  higher  elevations,  while  average  root  growth  is  the  opposite  [77]. 
The  early  rapid  growth  of  P.  brutia  continues  during  the  following 
years;  the  mean  annual  increment  can  be  over  10  m3/ha  in  appro¬ 
priately  spaced  plantations  between  the  ages  of  25  and  35  years,  as 
described  earlier  (site  class  I,  site  index  20  m  [35]. 

Calabrian  pine  is  the  main  forest  species  in  Cyprus  growing  from 
sea  level  up  to  the  1600  m  altitude  [79].  Its  forests  represent  about 
90%  of  the  Island’s  forest  area,  which  totals  at  1 73,1 82  ha  [5].  Most  of 
these  forests  were  classified  currently  unexploitable  after  a  recently 
carried  out  inventory  while  41,399  ha  were  classified  as  currently 
exploitable  [5].  Fellings  are  confined  in  the  latter.  The  annual  wood 
output  for  the  periods  1 988-1992, 1 998-2002  and  2003-2007  was 
56,100m3,  27,633  m3  and  14,384m3,  respectively  [5].  This  wood 
is  utilized  locally  for  the  production  of  sawn  timber,  chipboard, 
box-shooks,  fuelwood,  charcoal  and  other  minor  wood  products 
[79].  Sixteen  thousand  eight  hundred  thirty  five  hectares  of  the  P. 
brutia  forests  have  been  selected  and  preliminarily  classified  either 
as  National  Forests  Parks  to  provide  amenities  to  the  public,  or  as 
Nature  Reserves  for  nature  conservation  [5]. 

The  past  management  (silvicultural  system)  of  pine  forests  in 
Cyprus  proved  to  be  unsuitable,  resulting  in  serious  ecological, 
economic  and  social  constraints  [79].  This  may  put  in  danger 
the  survival  of  these  stands  in  the  long  term  [43].  The  degrada¬ 
tion  of  pine  forests  in  Cyprus  has  resulted  in  drastic  reduction  of 
annual  wood  output  from  56,100  m3  (for  the  period  1988-1992)  to 
14,384  m3  (for  the  period  2003-2007)  [5].  The  recognition  of  this 
situation  and  the  increasing  interest  in  using  these  stands  for  either 
direct  production,  or  indirect  production  (environmental  preserva¬ 
tion)  justifies  the  need  for  a  sustainable  management  of  Calabrian 
pine  stands.  Considering  the  high  environmental  and  silvicultural 
variability  of  P.  brutia  stands  [43],  it  is  necessary  to  typify  and  char¬ 
acterize  them  in  order  to  optimize  their  management. 

Estimating  forest  productivity  is  both  necessary  for  effective  for¬ 
est  management  and  useful  for  evaluating  basic  site  conditions  for 
ecological  field  studies.  Site  quality  is  therefore  influenced  by  fac¬ 
tors  such  as  available  light,  heat,  moisture,  and  nutrients,  along 
with  other  soil  characteristics  such  as  soil  depth  and  aeration  [29]. 
Although  it  would  be  best  to  directly  measure  and  predict  these 
factors,  they  require  precise  measurements  that  may  be  difficult 
to  extrapolate  across  scales.  Therefore,  indirect  methods  for  evalu¬ 
ating  site  quality  are  more  frequently  used  in  forest  management 
[56,67,72]. 

Site  index,  defined  as  dominant  height  at  some  fixed  base  age, 
is  one  of  the  most  commonly  used  indicators  of  site  productivity 
because  there  is  a  close  correlation  between  volume  and  dominant 
height  growth  [72].  Many  mathematical  functions  are  available 
to  model  dominant  height  growth.  Desirable  characteristics  for 
growth  functions  are:  polymorphism,  existence  of  inflection  point 
and  horizontal  asymptote,  logical  behavior,  right  theoretical  basis, 
base-age  invariance,  and  parsimony  [9,14,65].  These  requirements 
depend  on  both  the  construction  method  and  the  mathematical 


function  used  to  develop  the  curves  [65].  Among  other  meth¬ 
ods  for  site  index  curve  construction  [38],  the  algebraic  difference 
approach  presents  certain  advantages:  short  observation  periods 
can  be  effectively  used,  and  the  structure  of  equations  is  base-age 
invariant  [9,14,15,42]. 

Modeling  dominant  height  growth  for  Pinus  sp.  in  the  Mediter¬ 
ranean  region  has  received  some  attention.  Rio  and  Montero  [58] 
developed  site  index  curves  for  Pinus  sylvestris  in  Spain,  Calama 
et  al.  [67]  for  Pinus  pinea  in  Spain,  Bravo-Oviedo  et  al.  [2]  for  Pinus 
pinaster  in  Spain,  Hatzistathis  et  al.  [4]  and  Kitikidou  et  al.  [41  ]  for 
P.  brutia  in  Greece. 

Shater  et  al.  [85]  developed  dominant  height,  individual-tree 
diameter  increment,  tree  height,  and  self-thinning  site  curve  mod¬ 
els  forP.  brutia  stands  in  Syria.  In  this  work,  the  Chapman-Richards 
line  was  selected  to  show  the  development  of  dominant  height 
on  average  site.  Palahi  et  al.  [57]  developed  and  compared  two 
sets  of  models  (site  index  models  for  even-aged  and  uneven-aged 
stands)  which  enable  tree-level  simulation  of  the  development 
of  pure  and  mixed  stands  of  P.  brutia  in  North-eastern  Greece. 
Adamopoulos  et  al.  [73]  examined  the  impact  of  site  qualities  on 
P.  brutia  wood  characteristics  in  North-eastern  Greece,  while  Kav- 
gaci  et  al.  [6]  studied  the  post-fire  long-term  regeneration  of  the 
species  in  South-western  Turkey  through  a  floristic  survey.  Aert- 
sen  et  al.  [80]  compared  and  evaluated  five  different  modeling 
techniques  (multiple  linear  regression,  classification  and  regres¬ 
sion  trees,  boosted  regression  trees,  generalized  additive  models 
and  artificial  neural  networks)  to  model  site  index  of  homogeneous 
stands  of  important  tree  species  of  Turkey  (P.  brutia  among  others) 
by  using  environmental  variables  (soil,  vegetation  and  topograph¬ 
ical  variables)  as  predictors.  In  another  study,  Aertsen  et  al.  [81] 
compared  the  performance  of  the  same  five  modeling  techniques 
for  the  prediction  of  forest  site  index  (P.  brutia  stands  among  others) 
in  two  contrasting  ecoregions  (lowlands  in  Belgium  and  Mediter¬ 
ranean  mountains  in  Turkey).  Artificial  neural  network  models 
have  been  developed  to  estimate  volume  of  dominant  Calabrian 
pine  trees  [49,51].  Finally,  Diamantopoulou  et  al.  [50]  applied 
Principal  Components  Analysis  and  Hierarchical  Cluster  Analysis 
in  order  to  interpret  the  behavior  of  several  quantitative  vari¬ 
ables  of  reforestation  Calabrian  pine  trees  (dominant  height  among 
others). 

The  main  goal  of  this  study  was  to  develop  a  dominant  height 
growth  model  for  P.  brutia  growing  in  two  different  site  types  in 
central  Cyprus,  from  which  a  site  quality  model  would  be  created. 
To  achieve  this  objective,  the  variability  in  dominant  height  growth 
patterns  and  the  differences  among  the  sites  were  analyzed. 


2.  Materials  and  methods 

2.2.  Dataset 

The  study  was  conducted  in  central  Cyprus  (34°56N,  33°17E), 
in  middle  elevation  (300-750  m  a.s.l.)  forest  areas.  This  area  cov¬ 
ers  about  14,000  ha.  The  annual  rainfall  and  the  mean  yearly 
temperature  (for  1991-2000),  according  to  the  closest  meteoro¬ 
logical  station  (Kornos,  34°55N,  33°24E,  370  m  a.s.l.),  are  477.1  mm 
and  18.4°C  respectively.  The  geological  substrate  of  the  study 
area,  belongs  to  volcanic  sequence  diabase  dykes  with  pillow  lava 
screens  [30].  The  soil  texture  is  mainly  sandy  loam  to  sandy  clay 
loam  (Petrou  and  Milios  [62],  data  from  Forest  Department  of 
Cyprus). 

The  study  area  was  formerly  consisting  mainly  of  vine  and  olive 
yards.  The  majority  of  those  fields  were  abandoned  during  the  last 
century  and  nowadays  they  have  been  recolonized  mainly  by  P. 
brutia  trees  which  is  the  dominant  species  in  the  area.  A  few  olive 
yards  are  still  cultivated.  Moreover  three  built  up  areas  are  located 


Table  1 

Summary  statistics  of  sampled  trees  per  site. 
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Site 

No.  of  trees 

Mean  T 

SD  T 

T  min 

T  max 

Mean  H 

SDH 

H  min 

H  max 

1 

32 

27.69 

4.12 

24 

37 

13.06 

2.04 

9.90 

18.80 

2 

32 

90.09 

7.13 

76 

108 

8.93 

1.67 

5.70 

12.30 

H  (m)  is  total  tree  height  at  age  T  (years). 


in  the  study  area.  These  are  the  villages  of  Lythrodontas,  Mathiatis 
and  Kapedes. 

The  R  brutia  trees  form  stands  with  different  stand  structures 
and  variable  densities  in  the  wider  study  area.  In  some  cases  the 
R  brutia  trees  create  small  single  cohort  stands  [18].  The  R  brutia 
trees  in  those  stands  have  a  ground  cover  percentage  of  at  least  75%, 
forming  dense  stands.  Also  the  R  brutia  trees  can  be  classified  to 
dominant,  co-dominant  and  intermediate  according  to  their  crown 
class.  [18,22]. 

Particularly  the  abovementioned  single  cohort  stands  are 
located  in  two  different  site  types  as  regards  their  soil  depth  (site 
types  A  and  B).  The  soil  depth  in  the  plots  of  site  type  A  is  over 
100  cm,  while  the  corresponding  depth  in  the  plots  of  site  type  B  is 
about  30  cm.  The  soil  depth  was  measured  at  2  randomly  selected 
points  in  each  plot.  Thirty-two  plots  of  100  m2  (10  x  10  m)  were 
randomly  established  in  both  site  types  (16  in  site  type  A  and  16  in 
site  type  B).  From  each  plot  one  dominant  and  one  co-dominant  tree 
was  felled,  according  to  the  stratified  random  sampling  method, 
for  stem  analysis  purposes,  resulting  in  64  trees  (32  in  each  site 
type).  From  each  tree,  cross-sectional  discs  were  cut  and  removed 
from  the  stump  height  (^0m),  at  breast  height  (1.3  m)  and  every 
1  m,  up  to  the  bole.  The  last  disc  was  collected  from  the  3  cm-bole 
diameter.  In  each  cross-sectional  disc  the  number  of  annual  growth 
rings  was  counted  using  the  LINTAB  system  of  RinnTech  and  the 
program  TSAP-Win  [25].  In  stem  analysis  for  the  height  calcula¬ 
tion  the  modified  version  of  Carmeans’  algorithm  was  used  [40,82]. 
The  data  set  was  composed  of  3769  recordings  from  the  64  sam¬ 
pled  trees.  The  data  summary  of  the  sampled  trees  is  presented  in 
Table  1. 

Data  were  split  in  two  groups.  The  first  group  (80%  of  the  trees 
-  2905  recordings  from  51  trees)  was  used  as  fitting  data,  and  the 
second  group  (20%  of  the  trees  -  864  recordings  from  13  trees)  was 
used  as  validation  data.  This  method  for  testing  the  performance  of 
a  proposed  model  is  known  as  cross-validation  [21 ,29,52,59,60,70]. 

2.2.  Candidate  functions 

An  algebraic  difference  approach  has  been  used,  since  it  shows 
better  properties  and  performance  than  analogous  fixed-base-age 
equations  [13,15].  Four  models  were  selected  from  those  most  com¬ 
monly  used  in  forest  research  as  candidate  functions  to  model 
dominant  height  growth  (Table  2).  The  first  two  polymorphic  func¬ 
tions  derive  from  the  Bailey-Clutter  function  [66],  the  next  two 

Table  2 

Candidate  models  for  dominant  height  modeling. 


from  the  Korf  function  [78],  the  next  is  the  McDill-Amateis  func¬ 
tion  [55],  and  the  final  two  polymorphic  functions  derive  from  the 
Chapman-Richards  function  [24]. 

2.3.  Data  structure  and  model  fitting 

To  fit  an  algebraic  difference  equation  expressed  in  the  general 
form  of  H2  =/(Hi,  Ti,  T2),  different  data  structures  defined  in  Bor¬ 
ders  et  al.  [8]  can  be  used.  These  data  structures  are  relevant  in  any 
increment-based  modeling.  The  data  chosen  for  fitting  the  differ¬ 
ent  functions  included  all  the  possible  combinations  of  height-total 
age  pairs  for  a  tree  (all  possible  growth  intervals)  [8,28,39,67,74]. 

The  fittings  were  carried  out  using  the  nonlinear  procedure  on 
the  SPSS  software  [71  ].  The  Levenberg-Marquardt  iterative  method 
was  selected  because  it  is  the  most  useful  when  the  parameter 
estimates  are  highly  correlated  [84]. 

The  autocorrelation  derived  from  using  stem  analysis  data 
was  barred  by  applying  the  Goelz  and  Burk  [39]  correction.  First, 
each  function  is  fitted  following  ordinary  non-linear  least  squares 
regression  and  the  error  term  e^,  residual  from  estimating  using 
Hj,  is  expanded  following  an  autoregressive  process: 

eij  —  P£i-\,j  +  Y£i,j-\  +  £ij 

where:  p  =  autocorrelation  between  the  current  residual  and  the 
residual  from  estimating  H using  Hj  as  a  predictor  variable; 
y  =  relationship  between  the  current  residual  and  the  residual  form 
estimating  Hi  using  Hj_ t  as  a  predictor  variable;  Sjj  =  independent 
errors  with  mean  zero  and  constant  variance  rj2.  The  model  param¬ 
eters  are  then  obtained  by  fitting  the  expanded  function.  The 
autocorrelation  parameters  vary  the  weight  of  each  observation  by 
reducing  the  residual  proportional  to  a  previous  residual.  Neither 
of  the  autoregressive  parameters  p  nor  y  are  used  for  field  appli¬ 
cations  of  equations  because  the  errors  £i-\j  and  £ij-\  cannot  be 
observed  without  stem  analysis  [72]. 

2.4.  Model  selection  criteria 

A  three-step  procedure  was  used  to  evaluate  and  select  the  most 
appropriate  model,  which  included  qualitative  as  well  as  quantita¬ 
tive  examinations.  The  first  step  was  to  evaluate  the  model  fitting 
statistics  based  on  seven  model  performance  evaluation  criteria 


Model  Name 

1  Bailey-Clutter  [66] 

2  Korf  [78] 

3  McDill-Amateis  [55] 

4  Champan-Richards  [24] 


Height-age  equation 


H  =  ai  +CZ2F3 


H  =  a^a  2/ra3) 


H  = 


ai 

l+(a2/T«3) 


H  =  at[l  -e(-a2Da3] 


Free  parameter 

a-i 

a2 

a2 

a3 

a2 

a\ 

a2 


Difference  form 


H2  =  H,eb(T ?-r?3) 

\\2  _  ea1+(lnH1-a1)(T2/r1)a3 

H2=a1(H1/a,f>/7'2)°3 

Hi 


a\ 


ln  TV(1  /a2  )ln(72 /T! }  In  (aa  /fq  )ln 
al 


l-(l-(a1/W1  ))(T1/T2)a3 


u  u  (We^lY3 

u2  =  M  1  j^njry  J 


(1.1) 

(1.2) 

(2.1) 

(2.2) 

(3) 

(4.1) 

(4.2) 


Symbols:  H  is  total  tree  height  at  age  T;  H2  and  H 1  are  heights  at  age  T2  and  T\,  respectively;  a,  are  model  parameters  and  e  is  the  base  of  the  natural  logarithm  ln. 
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Table  3 

Model  performance  evaluation  criteria  (fitting  and  testing  procedures). 


Criterion  Symbol  Formula  Optimum  value 


Mean  residual 


Variance  ratio 


MRes 


VR 


E 


i= 1 


En  ~  ;  2 


0 


Residual  mean  of  squares 


RMS 


n-p 


i=l 


0 


Absolute  mean  residual 


Coefficient  of  determination-model  efficiency 


Linear  regression 


AMRes 

n  |  ^ 

|  A  —Hi 

0 

n 

1=1 

R2-MEf 

v— -  2 

>  M -Hi) 

i  Z — Ji=\ 

*  -  2 

1 

a,  /3 

Hi  =  a  +  /3Hi  +  Si 

a 

a  =  0,  /3  =  1 


Symbols :  H,  is  the  ith  estimated  value;  Hi  is  the  ith  observed  value;  n  is  the  number  of  observations;  p  is  the  number  of  model  parameters. 


described  by  Amaro  et  al.  [1  ]  (Table  3),  and  select  those  equations 
who  appeared  to  be  the  best. 

In  step  two,  the  characterization  of  the  model  error  was  ana¬ 
lyzed,  based  on  the  independent  data  set  testing  defined  in  Section 
2.1  (864  recordings  from  13  trees).  For  validation  purposes,  evalu¬ 
ation  criteria  applied  in  the  first  phase  were  also  calculated. 

Finally,  the  correctness  of  the  theoretical  and  biological  aspects 
of  the  models  was  assessed  in  step  3. 


Table  4 

Fit  statistics  for  each  model. 


Model 

MRes 

VR 

RMS 

AMRes 

MEf 

Linear  regression 

a 

P 

Ml 

4.0092 

0.0632 

25.5877 

4.0092 

-1.0147 

2.0971 

2.4969 

M2 

-0.0145 

0.3855 

7.6554 

0.0145 

0.3972 

-0.0954 

1.0153 

M3 

-0.0359 

0.3611 

7.7283 

0.0359 

0.3915 

-0.2602 

1.0421 

2.5.  Comparison  of  height  growth  models  among  sites 


Once  the  best  function  had  been  selected,  the  differences  in  the 
dominant  height  growth  models  for  the  different  sites  were  com¬ 
pared  using  both  the  full  and  the  reduced  models.  The  reduced 
model  corresponds  to  completely  different  sets  of  parameters  for 
different  sites  -  as  these  are  defined  from  the  selected  function  - 
while  the  full  model  corresponds  to  the  same  set  of  parameters  for 
all  the  sites  combined. 

Two  tests  for  detecting  simultaneous  homogeneity  among 
parameters  were  used:  the  Bates  and  Watts  non-linear  extra  sum 
of  squares  F  test  [74,75]  and  the  test  proposed  by  Lakkis  and  Jones, 
in  Khattree  and  Naik  [69],  to  compare  the  differences  in  site  index 
models  between  sites.  These  tests  are  frequently  applied  to  analyze 
differences  among  different  regions  [37,67,75]. 

Besides  the  full  and  reduced  models,  the  sum  of  squares  error 
(SS)  is  necessary  to  calculate  both  tests.  This  error  was  calculated 
as  follows: 


ss^EE(A-h,.)2 


n 


j= i 


where:  n  =  number  of  observations  for  each  tree;  m  =  total  number 
of  trees. 

The  F-test  is  applied  using  the  following  equation: 
(SSr-SSf)/(dfr-dff) 

C ssf)Mf ) 


The  L  statistic  used  in  the  Lakkis-Jones  test  is  defined  as: 


where:  SSf  and  SSr  =  error  sum  of  squares  for  full  and  reduced  model 
respectively;  and  m  =  total  number  of  trees.  If  homogeneity  exists 
among  the  model  vectors  of  parameters  jd,  the  distribution  of  the 
statistic  -21n(L)  converges  in  probability  to  a  Pearson  x2  distribu¬ 
tion,  with  v  degrees  of  freedom,  where  v  is  equal  to  the  difference 
between  the  number  of  parameters  estimated  in  the  full  model  and 
the  reduced  models. 


3.  Results 

Performance  results  of  all  models  are  shown  in  Table  4.  All  the 
parameter  estimates  for  all  the  functions  were  significant  at  an 
a  level  of  5%,  except  for  those  of  the  Champan-Richards  function 
(Model  4)  which  did  not  meet  the  convergence  criterion.  The  mean¬ 
ing  of  the  negative  value  of  MEf  for  the  Model  1  is  that,  the  mean 
value  of  observed  heights  gives  a  better  estimate  than  the  expected 
values  resulting  from  this  model  [1].  The  Korf  [78]  model  resulted 
to  be  the  best  model.  Differences  between  M2  and  M3  were  very 
small,  although  the  analysis  of  the  fit  statistics  revealed  that  the 
Korf  function  generally  resulted  in  slightly  lower  values  for  MRes, 
RMS  and  AMRes  as  well  as  higher  efficiency.  The  shape  of  the  curve 
of  the  Model  2  is  independent  of  the  autocorrelation  correction 
(Fig.  1  shows  both  fits  of  M2),  showing  practically  the  same  curve. 
Regression  coefficients  for  M2  are  calculated  as: 

H  =  10.924e(-5-876V(T°-648) 


where:  SSjandSSr  =  error  sum  of  squares  for  full  and  reduced  model  h  -  10  926e(_5-875)Ar°  647) 
respectively;  dfj  and  dfr  =  degrees  of  freedom  for  full  and  reduced 

model  respectively.  F  follows  an  F-distribution.  without  and  with  autocorrelation  correction,  respectively. 
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Table  5 

Testing  statistics  for  each  model. 


Model 

MRes 

VR 

RMS 

AMRes 

MEf 

Linear  regression 

a 

P 

Ml 

3.5235 

0.0892 

19.8586 

3.5235 

-0.8877 

2.2185 

2.1629 

M2 

-0.0015 

0.4404 

5.8193 

0.0015 

0.4468 

-0.0353 

1.0073 

M3 

-0.0092 

0.4393 

5.8113 

0.0092 

0.4476 

-0.0530 

1.0094 

Total  age  (years) 

Fig.  1.  Shapes  of  the  curves  of  Model  2  resulting  from  both  autocorrelation  correc¬ 
tion  and  not  autocorrelation  correction  fits. 

The  functions  M2  and  M3  were  again  those  which  per¬ 
formed  better  with  the  testing  data  (Table  5).  Once  more, 
Chapman-Richards  function  did  not  converged,  and  the  Models  2 
and  3  performed  similarly. 

The  shapes  of  the  curves  of  Model  2  seem  to  be  biologically 
reasonable,  which  prevent  unrealistic  height  predictions  when 
extrapolating  the  function  beyond  the  range  of  the  original  data. 
Model  2  was  selected,  based  in  the  good  results  in  the  performance 
criteria  defined  in  Table  3  with  both  fitting  and  testing  data. 

Assuming  the  suitability  of  the  Model  2,  it  would  be  interest¬ 
ing  to  analyze  the  dominant  height  growth  pattern  among  sites. 
Figs.  2  and  3  represent  the  Model  2  adjusted  for  each  site,  as  these 
are  defined  from  Model  2.  Site  I  shows  a  very  rapid  dominant  height 
growth  pattern,  especially  for  younger  ages. 

The  non-linear  extra  sum  of  squares  F  test  and  the  Lakkis-Jones 
test  revealed  that,  the  null  hypothesis  of  parameter  homogeneity 
was  acceptable  in  the  full  model  for  all  sites.  The  reduced  model 
for  site  I  did  not  converged  (Table  6).  On  the  other  hand,  the  fit 
statistics  obtained  with  full  and  reduced  models,  in  all  sites  and 

Table  6 

L  and  F  statistics. 


Fig.  2.  Site  index  curves  for  Pinus  brutia  in  central  Cyprus,  using  the  Korf  function 
(Model  2). 

sites  I  and  II  (productive  and  less  productive)  separately,  were  very 
dissimilar  (Table  7).  According  to  these  results,  the  total  model  (all 
sites  together)  can  be  selected.  In  our  attempt  to  apply  the  selected 

M2 

30  - 


0  20  40  60  80  100  120 


Age  (years) 

Fig.  3.  Site  index  curves  for  Pinus  brutia  in  central  Cyprus  using  the  Korf  function 
together  with  the  stem  analysis. 


M2 

df  full  model 

ssf 

df  reduced  model 

SSr 

Parameter  L 

Parameter  F 

All  sites 

2902 

22216.522 

2902 

22216.522 

1.000 

0.000 

Site  I 

Did  not  converged 

Site  II 

744 

1461.900 

2155 

3337.704 

0.0000* 

0.6766 

Significant  L-values  and  F-values  are  marked  with  *. 


Table  7 

Fit  statistics  of  the  full  and  reduced  models  for  the  two  sites. 


M2 

N 

Model 

MRes 

VR 

RMS 

AMRes 

MEf 

Linear  regression 

a 

P 

All  sites 

2905 

Full 

-0.0145 

0.3855 

7.6554 

0.0145 

0.3972 

-0.0954 

1.0153 

Reduced 

-0.0142 

0.3855 

7.6560 

0.0142 

0.3972 

-0.0943 

1.0151 

Site  I 

747 

Full 

0.0333 

0.9271 

1.9649 

0.0333 

0.9035 

0.1135 

0.9873 

Reduced 

Did  not  converged 

Site  II 

2158 

Full 

4.8979 

0.0514 

36.7628 

4.8979 

-0.8063 

0.3312 

4.1685 

Reduced 

1.2765 

1.2269 

23.3246 

1.2765 

-1.4401 

5.7317 

-0.2036 
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model  to  sites  I  and  II  separately,  the  reduced  model  for  site  I  did 
not  converged,  while  both  full  and  reduced  model  for  site  II  gave 
results  far  beyond  ideal. 

4.  Discussion  and  conclusion 

Accurate  estimates  of  forest  productivity  are  needed  for  sustain¬ 
able  forest  management  in  order  to  determine  annual  allowable  cut 
and  rotation  period,  and  to  make  tree  species  selection  decisions 
[11].  When  the  time  needed  to  produce  economically  feasible  and 
ecologically  sustainable  crops  is  unknown  and  existing  growth  and 
yield  information  is  limited,  developing  site  index  models  helps  to 
investigate  growth  performance  [11]. 

This  study  presents  a  site  index  model  for  the  P.  brutia  stands 
in  central  Cyprus.  The  algebraic  difference  approach  used  to  derive 
dynamic  height  functions  of  dominant  trees  assures  a  high  degree 
of  robustness  in  applications  [15].  According  to  fit  statistics,  the 
growth  function  of  Korf  [78]  was  chosen  to  explain  the  height 
growth  pattern  of  this  species.  The  statistic  criteria  with  testing 
data  showed  that  the  growth  function  of  Korf  was  still  the  best  pat¬ 
tern.  Testing  using  separate  intervals  (sites)  produced  worse  results 
that  working  with  all  possible  intervals  (sites),  so  it  seems  advisable 
to  analyze  all  data  together,  especially  in  the  testing  phase,  in  view 
of  the  fact  that  the  use  of  separate  intervals  is  not  very  common  in 
site  index  estimation  [23]. 

Comparing  the  results  of  the  selected  model  (M2)  among  sites, 
there  were  some  differences  in  error  measurements  and  model 
efficiency  values.  Similarly,  parameter  estimates  varied  consider¬ 
ably  among  sites.  This  could  be  due  to  little  data  available  for  an 
adequate  estimation  [16]  as  a  result  of  the  absence  of  old  stands, 
probably  because  of  their  previous  land  use  as  vine  and  olive  yards. 

An  effort  was  made  to  ensure  that  sampling  occurred  through¬ 
out  the  entire  range  of  Calabrian  pine  stands  in  the  study  area. 
Within  this  range,  it  seems  that  an  intermediate  site  is  missing 
(Fig.  3).  Although  Calabrian  pine  is  found  in  intermediate  sites  in 
the  study  area,  the  pine  trees  are  not  forming  so  dense  stands.  In 
this  sense,  the  general  population  of  pine  stands  was  sampled,  with 
canopy  cover  at  least  75%  of  the  ground,  as  mentioned  in  Section 
2.1.  Thus,  the  model  that  was  developed  provides  site  index  esti¬ 
mates  that  will  be  applicable  throughout  most  of  the  range  of  rather 
dense  pine  stands  in  the  study  area. 

In  spite  of  the  apparent  differences  observed  also  in  the  graph¬ 
ical  comparison  of  the  site  index  curves  (Figs.  2  and  3)  obtained  by 
fitting  the  selected  function  to  each  site,  the  statistical  tests  did  not 
reject  the  null  hypothesis  of  equality  of  height  growth  patterns. 

The  application  of  a  single  height  growth  model  for  all  the  stud¬ 
ied  area,  can  help  us  classify  the  quality  of  Calabrian  pine  stands  and 
decide  what  measures  can  be  taken  regarding  the  different  silvicul¬ 
tural  treatments  [37].  Nowadays,  due  to  the  state  of  degradation  in 
which  many  ofP.  brutia  stands  in  Cyprus  have  fallen  [79],  it  is  neces¬ 
sary  to  establish  priority  areas  where  resources  should  be  invested 
and  this  site  index  model  may  provide  a  key  tool  for  this  purpose. 
A  common  dominant  height  growth  model  could  also  simplify  the 
development  of  a  site  index  model  based  on  ecological  variables 
[3].  Although  it  is  sometimes  necessary  to  stratify  the  study  area 
when  looking  for  relationships  between  site  index  and  ecological 
variables,  the  analysis  is  facilitated  with  a  single  dominant  height 
growth  pattern  [29,31]. 

As  regards  comparison  of  P.  brutia  growth  patterns  in  the 
Mediterranean  region,  if  the  obtained  curves  in  this  work  are  com¬ 
pared  with  the  site  index  model  for  Greek  P.  brutia  stands  in  Thasos 
island  [41],  stands  in  central  Cyprus  show  higher  growth  rates  at 
younger  ages  and  a  larger  gap  between  the  best  and  the  worst  site 
qualities.  This  could  be  due  to  different  soil  depth  of  the  two  areas 
or  due  to  different  genetic  background.  Additionally,  the  large  sam¬ 
pled  area  of  this  study  and  the  good  results  obtained  in  different 


sites  (fitting  and  testing  data)  makes  this  site  index  model  a  good 
option  for  classifying  site  qualities  of  P.  brutia  stands  in  Cyprus. 
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